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Abstract 

A chaos control algorithm is developed to actively stabilize unstable periodic 
orbits of higher-dimensional systems. The method assumes knowledge of the 
model equations and a small number of experimentally accessible parameters. 
General conditions for controllability are discussed. The algorithm is applied 
to the Hamiltonian problem of point vortices inside a circular cylinder with 
applications to an experimental plasma system. 
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I. INTRODUCTION 


Hamiltonian models of 3-D collinear vortex filaments or 2-D point vortices inside closed 
geometries have for many years attracted interest due to its central role in understanding the 
evolution of vorticity in real flows 2-D point vortex dynamics is also relevant in a series 

of real experimental applications such as vortices in superfluid helium [7f, superconductors, 
and dynamics of non-neutral plasma filaments ||. It has also been shown recently that 
point-vortex models can provide a good formulation for developing control algorithms in 
realistic flows ||, since they are able to capture in certain limits most of the qualitative 
features of the vorticity dynamics. These results indicate that control algorithms based 
on these low-dimensional point vortex models can be successfully implemented for certain 
geometries even in fully-viscid Navier-Stokes equations 


In this Letter we will likewise attempt to develop a low-dimensional model for a fluid 
control algorithm, aimed at modifying a realistic flow of vorticity inside a bounded fluid 
domain. The model will be based on Hamiltonian point vortex dynamics, and make use of 
unstable periodic trajectories in the natural flow dynamics to modify the flow via the Ott, 
Grebogi, Yorke (OGY) JTT]] scheme. The principal result here is that we develop a general 
control algorithm which can be utilized in higher-dimensional systems (phase space D > 3), 
such as multi-vortex systems. We demonstrate the method’s effectiveness using numerical 
examples that model an experimental system of confined non-neutral plasma. 

To develop the control scheme, we first briefly review the formulation of point vortex 
dynamics. Let z = x + iy denote the complex coordinate of a single vortex. The effect of the 
boundary of the circular domain is simply accounted for by image charges at R 2 /z*. The 
Hamiltonian of a system of N identical point vortices each of circulation T inside a circle of 
radius R is then given by [[J 
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Each vortex moves by being passively advected in the velocity field created by the other 
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vortices and the image vortices, with the dynamical equations being obtained from the 
canonical Hamilton relations: 




dH 


k — 1, 


N, 


( 2 ) 


that is 


z u = 
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47T 
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j^k Zk z i j ~r 

Here, the hrst sum includes the velocity held created by the other vortices at the position of 

the fc-th vortex (inhnite self-interaction excluded) while the second sum gives the velocity 
held of the image vortices. The Hamilton relations (2) imply that each vortex contributes 
one degree-of-freedom and therefore two dimensions to the phase space. In addition, since 
the Hamiltonian does not contain time explicitly, the “energy” of the vortex system E = 
H(zi,...,Zn) is a constant of the motion. There exists another constant of motion for 
a circular domain resulting from the invariance of H under rotation, i.e. the “angular 
momentum”: 


N 

+ E 


R 2 - Z*Z k 


k — 1,..., N. 


(3) 


N 

£ = r£|a| 2 . (4) 

k 

Because of these two integrals of the motion, the dynamics of the iV-vortex system is re¬ 
stricted to a 2N — 2 dimensional manifold in the 2 N dimensional phase space. 

The Hamilton relations for point vortices imply the unusual situation that the phase 
space is simply a re-scaled version of the configuration space. Hence, integrability conditions 
which normally rely on symmetries of the Hamiltonian via Noether’s theorem, can partly 
be deduced by looking for symmetries in the physical domain M- In this view, one has 
global non-integrability (chaos) if N* > 1, where N* = (number of vortices) - (number 
of remaining integrals of motion) [1]. For the circular domain discussed here, N* > 1 for 
N > 3. In the following we will restrict our discussion to the three-vortex case, which is the 
smallest number of vortices inside a cylinder that leads to chaotic dynamics. 

Our central aim here is to achieve control and avoid chaos by stabilizing simple unstable 
periodic orbits of the three-vortex system, thereby taking advantage of the natural dynamics 
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of the flow. Although there are a large number of low-order periodic orbits of the generic 
vortex dynamics, we will restrict ourselves without loss of generality to those associated with 
symmetric configurations of vortices, with the center-of-vorticity in the origin. Examples 
of such configurations are shown in Fig. 1 for At = 3. We can simplify the discussion by 
rendering the equations-of-motion dimensionless with the substitution: 


z 

R' 


and t 


t 


2vr 

ER2' 


(5) 


Performing a further transformation into a reference frame co-rotating with the vortex sys¬ 
tem with uniform angular velocity u, again without loss of generality, Eq. (Q) will read 


as: 


Zu = —- 


E 


j^k Z k Z 3 j - -3 

The main advantage of this co-rotating frame is that some of the periodic orbits in the 6-D 
phase space reduce to a single fixed point (with proper choice of a;). 

The stability of such symmetric configurations has been first analyzed by Havelock W2 
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2loz1 


Z<Zk 


k = 1,2,3. 


( 6 ) 


In spite of the global non-integrability of the system, the symmetric orbits are stable for 
some parameter regimes, i.e. they are elliptic quasiperiodic orbits embedded in the four¬ 
dimensional manifold determined by the constants of motion (i and ((H). These orbits 
correspond to states where the vortices are away from the boundary of the cylinder, hence 
the effect of the “image” vortices is a weak perturbation, of the integrable three-vortex 
dynamics in free space. These regions are separated from the chaotic areas by KAM tori 
containing orbits of marginal stability. Many of these stable configurations have been found 
experimentally in the dynamics of magnetically confined non-neutral plasma columns [H, 


where the dynamics is very close to those of point vortices (13 . 


II. FORMULATION OF A HIGH-DIMENSIONAL CONTROL ALGORITHM 


Although many examples of chaos control algorithms have been developed in recent years, 
these have been restricted to a relatively low-dimensional phase space. Here, we develop a 
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control algorithm for higher-dimensional systems which assumes explicit knowledge of the 
dynamical equations and therefore the location, eigenvalues and eigenvectors of the fixed 
point or periodic orbit to be stabilized. This method can also certainly be reformulated in 
the language of experimentally measured time series (as in previous works [ |T4| -|T6| ), however 
this is not our immediate aim here. 

For the particular application we discuss in this Letter, we aim to modify the dynamics 
of the fully-viscid continuous fluid using a control model based on the point vortex dynam¬ 
ics. Hence, the model governing equations are known, and the high-dimensional control is 
straightforward to calculate. In the case of an experimental system, the method can be im¬ 
plemented in an identical fashion, using the Jacobian and the perturbation matrix obtained 
from experimental time-series m 


To begin, let r = A(r, S) represent the n = 2 TV-dimensional unperturbed dynamical 
system of Eq. (6). Here, in addition to Eq. @ we have introduced a set of m experimentally 
accessible system parameters, represented by the m-dimensional vector S. We assume the 
number of such parameters m is usually much smaller then the dimensionality of the phase 
space n, typically one. 

Suppose that r 0 is an unstable fixed point in the n-dimensional phase space, i.e., 
r 0 = A(r 0 ,$ 0 ) = 0, where So is a fixed set of parameters of interest. Allowing small per¬ 
turbations JS(£) in the parameters, then in a neighborhood of this fixed point r = r 0 + Jr 
the dynamics is described by the linearized equations 


Sr(t) = J(r 0 , S 0 ) Jr(£) + G(r 0 , S 0 ) JS(£) (7) 

where J = dA/dr is the standard Jacobian and the n x m matrix G = <9A/<9S describes 
the effect of small parameter perturbations on the system (i.e. the perturbation matrix). 

A stability analysis can be performed around the fixed point r 0 , by studying the proper¬ 
ties of the Jacobian, and typically can reveal a multitude of topologies. These topologies are 
defined by the set of eigenvalues of J, which can be real, complex or zero. The instability 
of the fixed point implies that there is at least one eigenvalue A with Re(X) > 0. Since 
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the system is Hamiltonian (TV J = 0) there is also at least one eigenvalue with negative 
real part. Although there is no guarantee of purely real eigenvalues and eigenvectors, with 
a proper transformation one can eliminate the imaginary part for at least one of the sta¬ 
ble eigenvalues (and the corresponding eigenvector) without changing the stability. In the 
following, we will suppose without loss of generality that there is at least one purely real 
negative eigenvalue A s < 0 with a corresponding normalized, stable real eigenvector e s . 

To achieve control of the vortex trajectories using small perturbations h<E>, we impose 
the condition that after a short time At the trajectory has approached the fixed point, 
i.e. | 5t( t + At) |<| 5r(t) \ . This can generally be accomplished in many different ways. 


One possible criteria is suggested by the low-dimensional chaos control of OGY |lT[] , where 
the new point 5r(t + At) is driven onto the stable manifold of the fixed point. Other 
possible choices are the self-locating (geometric control) method originally developed for 
low-dimensional chaos [l^|T^| or a method using a Newton algorithm developed for higher¬ 


dimensional systems |2(J. Here, we adopt the approach of OGY , i.e. we require that the 
trajectory lies on the stable manifold of the fixed point after a time At [9|: 


8r(t + At) = a 1 8r(t) \ e s 


( 8 ) 


where alpha is a small real number (0 < a < 1). Intuitively, this implies that after time 
At the trajectory lies on the stable manifold and simultaneously the distance from the fixed 
point has been decreased by a. Using a discretization of (0), this equation can be written 
explicitly as 


a | 8r(t) | e s ~ [1 + JAf] 5r(t) + G 8$>(t)At, (9) 

where for simplicity the notations J = J(r 0 , 3>o) an d G = G(r 0 , $ 0 ) have been introduced. 
This represents a system of n linear equations with m unknown perturbation parameters 
hT* that typically has no solution when m < n. This means that it is not possible to control 
such a system with one or a few control parameters in a one-step process, as described by 
Eq. (|). 
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Here, however we develop an alternative way to overcome this difficulty. Let us introduce 
the following notation for the RHS of Eq. (j9|): 

F(r, $) = [1 + JAt] dr(t) + G 6${t)At. (10) 

Then, instead of Eq. (|j) one can impose a weaker condition: namely, that after p steps the 
trajectory lies on the stable manifold, i.e. 

dr(t + pAt) = a \ dr(t) \ e s (11) 

or 

a |5r(t)|e s = F(F(...(F(r, ( 5#W),...),(5# (p - 1 )),5#^) (12) 

where the shorthand notation has been introduced for the perturbation at time t+kAt, 
i.e. $(t + kAt). Eq (p~2|) then can be written explicitly as 

M P-1 G h$ (1) + M p_2 G h$ (2) + ... + G h$ (p) = ^ (a | 5r(t) \ e s - M p dr{t)) (13) 

where M = 1 + JAt. This provides us with n linear equations and mp unknown 
perturbations. The number p is choosen as the smallest integer satisfying the condition 
mp > n. In the following, we will assume that there is a p such as mp = n. If n is not 
divisible by m, there are mp — n perturbation parameters that can be freely set to zero. 

We note that the size dr* of the linear neighborhood of the fixed point gives a natural 
limit for the allowed size of the perturbations, and we assume that none of the perturbation 
parameters <5$^ can be larger than a preselected value However, since Eq. (10) contains 
the yet unspecified control time parameter A t, the perturbations scale with this parameter 
as 1/At. This means that we should impose a limit on the product <5$ f,At rather than on 
alone. The upper bound of this product should naturally be set by the system’s physical 
limitations. For example, in the 3-vortex problem < hfcAt is a length that is naturally limited 
by the typical size of the linearized region 5r*, i.e. h<3>*At = dr*. 
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It is possible to derive general conditions for the controllability of the system. To have a 
solution for the p-step control process (13), the square matrix C formed from the p matrices 
of size n x m, i.e. 

C = {m^GiM^Gj.-.jG} (14) 

has to have nonvanishing determinant. In two-dimensions, the system is always controllable 
provided the vector G is not collinear with e s . Intuitively, this means that small perturba¬ 
tions of the system do not move the fixed point exclusively along the stable direction. In 
higher-dimensions, such a condition is not strong enough; the perturbation may access only 
a very limited subspace, or may have a certain symmetry that leads to un-controllability, 
even if the two-dimensional analog condition is satisfied. Instead 

det C ^0 (15) 

incorporates all these additional conditions, and can be viewed as the basic criteria of con¬ 
trollability. This relation is also a useful tool to select appropriate perturbations for each 
particular control problem, and to find the minimum number of perturbations that can suc¬ 
cessfully control the system. The control parameter selection is probably the most important 
step in any such problem, especially when the Jacobian is highly-symmetric as in the case of 
symmetric periodic orbits of identical vortices. For such orbits it is straightforward to show 
that most of the symmetric perturbations lead to un-controllability. An example of this is 
the simple perturbation of “squeezing” the boundary, i.e. making it slightly elliptic with a 
small eccentricity. Here, the equations of motion can be simply deduced by using a simple 
conformal mapping f(z) = z + ez 3 that maps the unit circle into an ellipse, with semi-axis 
1 + e and 1 — e, respectively (for e<l). A careful analysis of the collinear configurations of 
Fig. 1(a) shows that the determinant of the corresponding controllability matrix C vanishes, 
due to the symmetry of the Jacobian J and the perturbation matrix G with respect to the 
non-central vortices. This means that it is not possible to control this particular configura¬ 
tion by using cylinder squeezing. In fact such behaviour is expected to be quite general in 



systems of interacting identical subsytems that respond in the very same way to external 
perturbations. 


III. NUMERICAL EXPERIMENTS 

To illustrate the above control algorithm, we actively stabilized the two simplest unstable 
periodic orbits of the three vortex system inside a cylinder, i.e. the symmetric states shown 
in Fig. l(a,b). For the collinear state of Fig. 1(a), the angular velocity of the configuration is 
a j— (3+a 4 )/(l — a 4 )/(4a 2 ), while for the triangular state of Fig. 1(b), it iso; = (l±2a 6 )/(l — 
a 6 )/(2a 2 ). Here a is the radial distance of the non-central vortices from the center of the 
cylinder. To simplify the equations of motion, the dynamics is viewed in a reference frame co¬ 
rotating with the vortices with this constant angular velocity u. In this frame, the periodic 
orbits become single fixed points. Analytic calculation of the Jacobian for the collinear state 
at a = 0.5 reveals the existence of four real eigenvalues, Ai i2 = ±4.7565 and A 3i4 = 0, and 
two purely complex ones, A 5 6 = 1.6027b The triangular state at a = 0.7 leads to two purely 
real eigenvalues \\.2 = ±0.0001165 and four complex ones, A 3 , 4 ) 5,6 = ±0.6168 ± 0.1670b 
The stable eigenclirection e s is then obtained in each case as the eigenvector associated with 
the single negative real eigenvalue. One can observe that although the two topologies are 
completely different, the controller requires knowledge of only one of the possibly many 
stable directions. 

To control the vortex dynamics we require non-symmetric and non-homogeneous per- 
turbatios. We have introduced simple perturbations applied only at the cylinder surface: 
namely, a uniform distribution of sources and sinks with variable strength s = So cos (kd). 
Although such a perturbation may be difficult to implement in a real fluid experiment, it is a 
real possibility in the plasma analog of the vortex dynamics [|J, where similar perturbations 
can be generated by external electric fields. For k = 1,2, and 3, this leads to a uniform, a 
quadrupole, and a sextapole held, respectively, as shown in Fig. 2(a-c). Since all these fields 
are highly symmetric, none of them alone can effectively control the vortex dynamics. How- 


9 


ever, their linear combinations, resulting in strongly asymmetric and inhomogenous fields, 
proved to be a successful choice that satisfies the controllability conditions (jT5p. Here, we 
have introduceed three linearly independent external fields as shown in Fig. 2(d-f). Using 
these, the equations of motion can be written: 


U = 


E 


i+k Zk z 0 j 

2 


+ E 


2 LUzt 


Z jZ k 


+<f>i (1 + Zk + z k ) + $2 (—1 + Zk — z k ) + $3 (1 + Zk — z k ) k — 1 , 2 , 3. 


(16) 


We previously showed that there are two constants of motion in the unperturbed sys¬ 
tem, the energy and the angular momentum. Due to the presence of the time-dependent 
perturbations that do not preserve these first integrals, however, the effective dynamics is 
no longer restricted to a 4-D manifold, but rather can occupy in the full 6 -D phase space. 
Note that if the dynamics initially lies in the same energy and angular momentum manifold 
that contains our fixed point, ergodicity ensures that the trajectory sooner or later reaches 
a close neighborhood of the desired fixed point, and then the control algorithm is switched 
on. 


If this is not the case, then in general this procedure can be supplemented with a targeting 
algorithm that first drives the system close to the proper energy and angular momentum val¬ 
ues, then takes advantage of the infinite number of unstable periodic orbits in that manifold 
to reach the neighborhood of the desired fixed point in a short time |H[] . The development 
of such a targeting algorithm in higher dimensions is clearly an important step that deserves 
further study. 

Since this is not the primary emphasis of this Letter, in our numerical simulations we 
have started the system close to the unstable fixed point. To demonstrate this Fig. 3(a-b) 
shows the time evolution of a single coordinate pair (xi,y i) of the vortex system in the 
collinear case, at first without control (up to mark A). Then, after the dynamics reaches 
a small neighborhood of the desired fixed point (mark B), the controller is switched on, 
and the trajectory is stabilized. Fig. 3(c-e) also shows the time evolution of the required 
perturbations. Clearly this demonstrates that the algorithm can effectively control the 
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dynamics with only tiny perturbations applied on the boundary. For better visualization, 
Fig. 4 shows a three-dimensional projection of the phase space trajectory with and without 
control. 


As previously mentioned, other control criteria can be used instead of the OGY method. 
For example, one can consider the high-dimensional control of Xu and Bishop, based on 
the Newton root finding algorithm 0. This method does not assume knowledge of the 
stable direction, and it can be useful in cases where all eigenvalues are complex (with at 
least two having non-vanishing real part). Then, instead of condition (11) one can use the 
usual Newton root finding algorithm 


r (t + pAt) = r (t) - J (r 0 ) A(r 0 + hr). 


(17) 


However, this method fails when the Jacobian is not invertible. (We note here that the 
transition from Eq. (5) to Eq. (6) in Ref. is mathematically incorrect, i.e. the solution 
(6) does not fulfill Eq. (5). Therefore the weaker condition ( |T7| ) shown above should be 
used.) 

In a realistic Navier-Stokes simulation our algorithm is effective on time scales signifi¬ 
cantly shorter than the viscous one (Re 1 / 2 ). This is attractive since this time scale can be 
rather large in some plasma systems, where the viscosity is small and there are no boundary- 
layer effects due to the free-slip conditions on the boundary. Even in these circumstances, 
it should be noted that there is an additional non-viscous effect that is not present in the 
Hamiltonian model, the vortex merger due to the finite size of the vortices. Therefore ef¬ 
fective control can be reached only with concentrated vortices that are far from each other, 
typically with 2 p/d ^ 1.8, where d is the distance between vortices and p the vortex radius 
computed according to Ref. |22j]. Although this condition can be easily fulfilled while the 
controller is working, during the targeting algorithm or full chaotic dynamics the vortices 
may easily come close to each other. In this context, the principal usefulness of such a 
control scheme will be to prevent the vortex merger on short time scales from symmetric 
unstable initial conditions. 
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In conclusion, we have demonstrated that the control algorithm derived in this Letter 
can effectively control higher-dimensional systems. We have also derived a controllability 
condition that does not depend on the particular control method used, since it contains only 
the Jacobian and the perturbation matrix as input parameters. In this sense, the result is 
quite general, and can be used to decide the type and number of perturbations needed to 
control a higher-dimensional system. 

For the particular numerical example presented here, i.e. three-vortex dynamics inside a 
cylinder, the controller was designed in such a way that it could be implemented in a mag¬ 
netically confined plasma experiment. At present, we also plan to implement the control 
algorithm in a realistic viscous fluid framework, and hence an effective targeting algorithm 
has to be developed. In addition, in the present paper the controller was formulated using 
exactly known vortex coordinates. However, the control scheme itself can be entierly refor¬ 
mulated in a phase space reconstructed purely from wall signal (e.g. boundary pressure or 
voltage) measurements, thereby providing an algorithm which can be directly implemented 
in an experimental framework. 

After the acceptance of this paper for publication it was brought to our attention a 
different formulation of the control algorithm that yields to a similar controlability matrix 
23|. 
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FIGURES 

FIG. 1. Symmetric vortex configurations leading to periodic orbits. Three-vortex configura¬ 
tions: (a) collinear with one vortex in the center; (b) equilateral triangle. This configuration leads 
to stable dynamics for a < 0.567, and unstable one for a > 0.567. Each configuration is rotating 
with an angular velocity lu that depends on the particular configuration. 

FIG. 2. Perturbation fields used to control the vortex dynamics. Figs, (a-c) show the velocity 
fields generated by a distribution of sources and sinks on the boundary of strength cos ( k6 ) for 
k = 1,2,3, respectively. Figs, (d-f) displays the velocity fields resulting from different linear 
combinations of the fields (a-c): (d) z* = 4>i (1 + z + z 2 ), (e) i* = < 4>2 (— 1 + z — z 2 ), and (c) 
z* = d >3 (1 + z — z 2 ). The strength of these fields &k,k = 1,2,3 are computed by the control 
algorithm. 

FIG. 3. (a) Vortex trajectories in the lab frame and the applied perturbations. Figs, (a-b) 
show the coordinates of one of the vortices, with the collinear configuration started close to the 
fixed point. Control is off until time t = 20 (mark A). At this point the controller is activated 
and control begins at mark B. Without a targeting algorithm, the time between A and B can be 
very long. From instant B, the controller applies the necessary perturbations at intermittent times 
At = 0.1 to keep the trajectory on the unstable periodic orbit. The parameter a in Eq. (13) was 
0.1 . Figs, (c-e) show the applied perturbations. After a short transient time the perturbations 
are quite small, invisible on this scale. Fig. (e), inset, also displays the perturbation $3 magnified 
1000 times between 30 < t < 40, for visualization. 

FIG. 4. Three-dimensional projection (aq, y\. r \2 =| £2 — z\ |) of the phase space trajectory 
in the lab frame, (a) Uncontrolled trajectory of the collinear configuration, (b) Trajectory after 
controller is switched on. 
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